Method for Denoising the Vibration Signal of Rotating Machinery through VMD and MODWPT

The vibration signals from rotating machinery are constantly mixed with other noises during the acquisition process, which has a negative impact on the accuracy of signal feature extraction. For vibration signals from rotating machinery, the conventional linear filtering-based denoising method is ineffective. To address this issue, this paper suggests an enhanced signal denoising method based on maximum overlap discrete wavelet packet transform (MODWPT) and variational mode decomposition (VMD). VMD decomposes the vibration signal of rotating machinery to produce a set of intrinsic mode functions (IMFs). By computing the composite weighted entropy (CWE), the phantom IMF component is then removed. In the end, the sensitive component is obtained by computing the value of the degree of difference (DID) after the high-frequency noise component has been decomposed through MODWPT. The denoised signal reconstructs the signal’s intrinsic characteristics as well as the denoised high-frequency IMF component. This technique was used to analyze the simulated and real-world signals of gear faults and it was compared to wavelet threshold denoising (WTD), empirical mode decomposition reconstruction denoising (EMD-RD), and ensemble empirical mode decomposition wavelet threshold denoising (EEMD-WTD). The outcomes demonstrate that this method can accurately extract the signal feature information while filtering out the noise components in the signal.


Introduction
Rotating machines are important types of equipment used in a variety of industrial applications. An unexpected fault in these machines could lead to significant economic losses and casualties. The vibration signal of rotating machinery carries information about its running state. Therefore, detection and analysis of the vibration signals of rotating machinery is usually the basis for understanding the service life and operating status of the equipment. It is an important basis used for preventive maintenance and fault diagnosis [1]. However, due to the interference of operating equipment and the field environment, various kinds of noise are inevitably introduced in the process of signal monitoring and acquisition. In order to ensure the authenticity of the measured signal and the effectiveness of subsequent signal feature extraction, it is particularly important to denoise the measured signal.
Currently, there are many denoising methods used in the analysis of vibration signal of rotating machinery, including the Fourier filter, wavelet transform (WT), and fast independent component analysis (FastICA) [2][3][4]. While these methods have achieved some success, they often encounter problems. Traditional denoising methods are mainly based on linear filtering [5], but vibration signals of rotating machinery are often non-linear and non-stationary due to environmental noise and varying operation states [6]. As a result, traditional methods are not ideal for handling such signals. WT also has some unavoidable defects, such as boundary distortion, energy leakage, and non-adaptiveness [7]. For the Table 1. The characteristics of common signal decomposition methods.

Advantages Disadvantages
WT Multi-resolution analysis. Non-adaptive.
MODWPT Translation invariance of wavelet coefficients and scale coefficients.
Unreasonable decomposition effect of fast-changing signals.
Based on VMD and MODWPT, this paper proposes a novel signal denoising method. This method is used to denoise simulated signals and measured rotating machinery fault signals. The main contributions of this paper are summarized as follows: (1) The number N of effective center frequencies and minimum value of multi-scale dispersion entropy (MDE) are proposed to determine the key parameters K and α of VMD; (2) The comprehensive weighted entropy (CWE) of IMF components is used to select the sensitive IMFs of VMD, containing the signal feature information. The highfrequency noise IMFs are decomposed through MODWPT and the noise components are filtered out by calculating the difference degree (DID). The non-stationary rotating machinery signals are denoised by reconstructing the sensitive components of VMD and MODWPT; (3) The effectiveness and performance of the proposed method is verified by analyzing the simulated signal and real experimental gear fault signal. These results are compared with wavelet threshold denoising (WTD), EMD reconstruction denoising (EMD-RD), and EEMD wavelet threshold denoising (EEMD-WTD) to evaluate the performance.
The rest of this paper is organized as follows. Section 2 provides a brief overview of the fundamental theories. Section 3 describes the calculation process for the proposed denoising method. Section 4 provides a simulated signal example to demonstrate the reliability and effectiveness of the proposed method. Section 5 discusses the use of this method in fault detection and denoising of vibration signals from rotating machinery. Section 6 contains the paper's main conclusions.

VMD Method
The VMD method can decompose a signal f into a series of modal functions u 1 , u 2 , . . ., u k according to the preset scale parameter K. The resulting constrained variational problem is the following: where {u k }: = {u 1 , . . ., u K } and {ω k }: = {ω 1 , . . ., ω K } are shorthand notations for the set of all modes and their center frequencies, respectively. ∂ t is the partial derivative of the time t for the function, δ(t) is the unit pulse function, and * represents the convolution operation.
We introduce the augmented Lagrangian ζ as given below, and the constrained variational problems are transformed into unconstrained variational problems: The solution to the original minimization problem is now found as the saddle point of the augmented Lagrangian in a sequence of iterative sub-optimizations called the alternate direction method of multipliers. The multipliers include u n+1 k , ω n+1 k , and λ n+1 , where u n+1 k is the modal function at the (n + 1)th cycle, ω n+1 k is the center frequency of the power spectrum of the current modal function, and λ n+1 is the multiplication operator at the (n + 1)th cycle.
Then, the modal component u k and the center frequency ω k are derived: whereû n+1 k ,f , andλ n+1 represent the Fourier transform corresponding to u n+1 k , f , and λ n+1 , respectively. For convergence accuracy e > 0, when Equation (5) is satisfied, the decomposition stops. Thus, the final modal componentû k and its center frequency ω k are given as

VMD Key Parameter Selection Method
The preset scale parameter K and penalty factor α directly affect the accuracy of the VMD decomposition results. Although the current intelligent search algorithm can optimize the parameter values given above, such methods take a long time and it is difficult to achieve the purpose of actual detection. Therefore, this paper proposes an efficient and simple method of VMD decomposition key parameter selection.

Selection of K
The center frequency of each IMF component obtained through VMD is distributed from a low to a higher level. If the optimal value of K is obtained, this means that the center frequency distribution between adjacent IMF components is reasonable and the final result value will not be similar or mixed. Therefore, researchers currently use the center frequency observation method to determine the optimal value of K [29,30]. However, these methods lack quantitative judgment standards and universality in the actual signal analysis process.
From Equation (5), it can be seen that VMD uses the frequency domain change of the IMF component during two consecutive cycles as the iterative constraint condition to determine when to stop the decomposition. To solve this problem, this paper proposes an algorithm to select K according to the number N of effective center frequencies by taking the change in the center frequency value between two adjacent decomposition processes of the same order modal component as the judgment condition. Figure 1 shows the algorithm flow chart. The specific process of the algorithm is as follows: Step 1: Initialize the value of K, let K = 2.
Step 2: Decomposed the signal through the VMD method and obtain K order IMF components and the center frequency of each IMF component ω K,i (i = 1, 2, . . ., K).
Step 4: According to Equation (6), calculate the accuracy of the center frequency of each IMF component under the same order and different decomposition times ε K,k (k = 1, 2, . . ., K): where ω K+1,k is the center frequency of the (K + 1)th order IMF component obtained through VMD and ω K,k is the center frequency of the K-th order IMF component obtained through VMD.
Step 5: Determine the size of the judgment accuracy ε K,k and the accuracy threshold θ. If the value of θ is too small, it is difficult to effectively distinguish the frequency domain characteristics between components. If the value of θ is too large, it is also difficult to accurately reflect the frequency domain difference between each component because the vibration signal of rotating machinery is often composed of multiple frequency components. Therefore, within a reasonable range, the value of θ selected should be as large as possible. Based on [31] and a large number of tests, the value of θ is set to 0.15. If ε K,k ≤ θ, it is considered as the effective center frequency. If ε K,k > θ, it is regarded as the invalid center frequency. Let N K+1 and N K be the number of consecutive effective center frequencies obtained from the IMF1 component, whereas K + 1 and K are used to determine the number of mode functions.
Step 7: If N K+1 ≤ N K , the number of consecutive effective center frequencies no longer increases, so it is considered that there is suspected over-decomposition, and thus the provisional optimal value is K.
Step 8: Consider that the current provisional optimal value is K. If the signal to be decomposed contains multiple and different frequency components, there is still the possibility of under-decomposition. To avoid the above problems, let K = K + 2 and repeat step 3-step 5. If N K+2 ≤ N K+1 , the optimal preset scale is K; if N K+2 > N K+1 , then repeat step 3 to step 7 and obtain the second provisional value. If N K+3 ≤ N K+2 , the second tentative optimal value is K + 2, and the parameter corresponding to max{N K , N K+2 } is the optimal value.

Selection of α
α is another important parameter to be set during the process of VMD decomposition, which determines the bandwidth of IMFs. The larger the α value, the smaller the bandwidth of the IMF component. As the VMD algorithm shows good noise robustness, after the vibration signal of rotating machinery is decomposed, the components of background interference and environmental noise in the signal should be filtered out. Additionally, the reconstructed signal should contain more feature information and demonstrate strong regularity and self-similarity [6,21].
Multi-scale dispersion entropy (MDE) is an important factor in evaluating the complexity of non-stationary signals [32,33]. MDE is based on dispersion entropy (DE); compared with multi-scale sample entropy (MSE) and scale arrangement entropy (MAE), MDE effectively simplifies the calculation process, enhances the accuracy of feature extraction, is less affected by abrupt signals, and has higher algorithm stability. Based on the analysis above, if VMD is optimized with an optimal α, it will enhance the regularity and self-similarity of a rotating machinery's fault signal. Consequently, the MDE value of the reconstructed signal will be minimized.
For a signal {x(t), t = 1, 2, . . ., T}, a compound coarse-grained signal is used and the k-th coarse-grained sequence under the set scale factor τ is x τ k . The coarse-grained signal can be obtained as k,2 , · · · . Calculate DE with each scale factor τ of the coarse-grained signal, and MDE can be formed as where m is the dimension of the embedded vector, c is the number of categories, and d is the time delay. After selecting K, calculate the MDE of the reconstructed signal obtained through VMD of the penalty factor α in different value ranges. The selection steps for parameter α are outlined as follows: Step 1: Initialize the value of α; let α = 100.
Step 2: Calculate the MDE of the reconstructed signal obtained through VMD and obtain MDE1.
Step 3: To ensure the bandwidth of IMF components, let α = 200, calculate the MDE of the reconstructed signal obtained by VMD, and obtain MDE2.
Step 5: Select α corresponding to min{MDEi} as the optimal value.

Illusive IMF Component Selection Algorithm Based on Comprehensive Weighted Entropy
At present, there are few methods for identifying illusive IMF components. Existing methods often rely solely on signal characteristics in the time or frequency domain, making it difficult to comprehensively and accurately judge. Therefore, this paper proposes an illusive IMF identification method based on CWE.
In rotating machinery equipment, when a fault occurs, the energy in the fault signal frequency band will significantly differ, and energy entropy (EE) can effectively describe the change in signal energy with frequency distribution [34]. Power spectrum entropy (PSE) is a nonlinear characteristic quantity that characterizes the complexity of the signal. It can also characterize the distribution of the vibration spectrum type of the signal in the frequency domain [35]. The time-amplitude-frequency product can effectively reflect the change relationship between time, energy, and frequency and highlight the local characteristics of the signal, and has good time-frequency resolution. Time-amplitude-frequency product entropy (TPE) fully utilizes the advantages of entropy in signal information evaluation while avoiding the influence of similarity between fault features. It can more effectively depict the internal characteristics of the signal and has more accurate signal feature extraction ability [36]. The specific process of the illusive IMF component selection method based on CWE is as follows: Step 1: For each IMF component {u i (t), t = 1, 2, . . ., T}, calculate the energy e(u i ), power spectrum p(u i ), and time-amplitude-frequency product q(u i ).
where X(u i ) is the discrete Fourier transform of each sample in signal component u i (t) and S j (j = 1, 2, . . ., m) is the energy of each time-frequency plane of the Hilbert spectrum (the entire time-frequency plane is divided into m blocks, on average). Step 2: Normalize the energy, power spectrum, and time-amplitude-frequency product of all IMF components.
where K is the number of IMF components decomposed through VMD.
Step 3: Calculate the entropy increment of each IMF component.
Step 4: Obtain the CWE value of each IMF component according to Equation (12).
Step 5: Calculate the principal component factor λ i based on CWE.
Step 7: Calculate the λ i difference of two adjacent components.
Step 8: Find the index m of the maximum difference. Then, the first m components of sequence {λ } contain the main characteristic information of the signal and represent the sensitive IMF component. The IMF component corresponding to λ i = 0 is the noise or illusive component, which should be eliminated.
The proposed method comprehensively considers the relationship between each IMF component and the decomposed signal in the time domain, frequency domain, and timefrequency domain. This method can enhance the main components of the signal, weaken the interference components irrelevant to the signal characteristics, and ensure the purity of the signal characteristics. Moreover, it avoids the need to propose a discrimination threshold, and the accuracy of the calculation results is more scientific and reliable. Figure 2 shows the algorithm flow chart.

MODWPT Method
Generally, the IMF component that characterizes the high-frequency noise component contains not only noise but also useful information about the signal [37]. If these high-frequency components are directly filtered as noise interference, key information in the signal may be lost, reducing the accuracy of the signal feature extraction. In order to enhance the signal-to-noise ratio and the denoising effect, the MODWPT method is used to decompose the IMF component that characterizes high-frequency noise.
The decomposition coefficient of MODWPT can be represented by Wj,n = {Wj,n,t, t = 0,…, N − 1}, where j is the number of decomposition layers, n can be regarded as a frequency index that changes with j, and the decomposition coefficient of MOWDPT can be calculated as follows: When a rotating machine malfunctions, usually only a certain component malfunctions, such as the rotor, gear, or bearing. Therefore, only a small part of the vibration information can be used to indicate the failure of the rotating machinery after using MOT-WPT to decompose the high-frequency noise components containing some information on the fault characteristics. The fault characteristics are concentrated only on a certain component. At the same time, there is a certain correlation between the fault characteristics and the signal characteristic information under normal conditions. The noise interfer-

MODWPT Method
Generally, the IMF component that characterizes the high-frequency noise component contains not only noise but also useful information about the signal [37]. If these highfrequency components are directly filtered as noise interference, key information in the signal may be lost, reducing the accuracy of the signal feature extraction. In order to enhance the signal-to-noise ratio and the denoising effect, the MODWPT method is used to decompose the IMF component that characterizes high-frequency noise.
The decomposition coefficient of MODWPT can be represented by W j,n = {W j,n,t , t = 0,. . ., N − 1}, where j is the number of decomposition layers, n can be regarded as a frequency index that changes with j, and the decomposition coefficient of MOWDPT can be calculated as follows: where if n mod4 = 0 or 3, r n,l = { g l }, if n mod4 = 1 or 2, and r n,l = { h l }. Among them, g l and h l represent the scale filter and wavelet filter of MODWT, respectively. When a rotating machine malfunctions, usually only a certain component malfunctions, such as the rotor, gear, or bearing. Therefore, only a small part of the vibration information can be used to indicate the failure of the rotating machinery after using MOTWPT to decompose the high-frequency noise components containing some information on the fault characteristics. The fault characteristics are concentrated only on a certain component. At the same time, there is a certain correlation between the fault characteristics and the signal characteristic information under normal conditions. The noise interference component contains a few normal signal characteristics; hence, there is a large degree of information difference between the noise interference component decomposed through the MODWPT and the normal signal.
The model formula of the DID [38] can be expressed as follows: where µ n (i), δ n (i), µ un (i), and δ un (i) are the mean and standard deviations of the normal and fault signals, and x i and y i are normal and fault signals. By combining the DID technique and the MODWPT method, the noise interference components in the high-frequency IMF component can be effectively filtered out.

The Signal Denoising Method Based on VMD and MODWPT
The vibration signal of rotating machinery is always mixed with different noises during the acquisition process, which affects the accuracy of the signal feature extraction. The steps of the denoising method based on VMD and MODWPT are as follows: Step 1: Use the VMD method to decompose the vibration signal of rotating machine, and determine the optimal value of K based on the successive number of effective center frequencies.
Step 2: After selecting K, calculate the MDE of the reconstructed signal obtained through VMD of the penalty factor α in different value ranges. Select α corresponding to min{MDE} as the optimal value. Step 3: Calculate the principal component factor λ i based on the CWE of each IMF component, determine the illusive IMF components, and obtain the preliminary denoising signal x (t).
Step 4: Based on the center frequency of the illusive IMF components, decompose the IMF components (λ i = 0), which represent the high-frequency illusive component, through MODWPT. Calculate the value of the DID between the component obtained through MODWPT and the normal signal, filter out the noise components, and determine the components containing the characteristic information of the signal to further improve the denoising effect.
Step 5: Reconstruct the IMF component denoised through MOWDPT and the IMF components of each order characterizing the signal's characteristics to form a denoised signal x (t).
The algorithm flow chart is shown in Figure 3.

Simulation Analysis
Usually, the vibration signal of rotating machinery can be simulated using the superposition of the amplitude modulation, frequency modulation and the Gaussian white noise [39]. Based on the research findings of Sun et al. [39], the following analogue signals are established to more realistically simulate the multi-component vibration signals of rotating machinery:

Simulation Analysis
Usually, the vibration signal of rotating machinery can be simulated using the superposition of the amplitude modulation, frequency modulation and the Gaussian white noise [39]. Based on the research findings of Sun et al. [39], the following analogue signals are established to more realistically simulate the multi-component vibration signals of rotating machinery: where n(t) is the Gaussian noise signal. It is obtained by using the 'randn' function in MATLAB software (version 2016a) and its SNR is −9.5 dB.
The sampling frequency is 4096 Hz and the sampling time is 1 s. The time domain waveforms of the simulation signals x(t) and z(t) are shown in Figures 4 and 5.  From Figure 4, it can be seen that the simulated signal contains three main center frequencies: 170 Hz, 540 Hz, and 910 Hz. When the rotating machinery signal is weak or covered by strong background noise, the signal characteristics (impact property and center frequency) may not appear apparently. Comparing Figure 4 with Figure 5, it can be seen that the periodic impact components in the simulated signal are completely submerged by noise, the center frequency 910 Hz is difficult to distinguish, and the highfrequency part is severely affected by noise, with many interference spectral lines appear in this part.
In order to denoise the simulated signal, the VMD-MODWPT method is used. Firstly, the simulated signal z(t) is decomposed through VMD. The center frequency of each IMF component of z(t) is obtained under different K values. In order to obtain the optimal value of K, the judgment accuracy of each IMF component corresponding to different K values is calculated according to the method proposed in this paper. The calculation results are shown in Table 2. Table 2 shows that N4 = 4 > N5 = 1, so K = 4 can be tentatively set as the optimal value and N6 = N5 = 1. Therefore, the optimal number of mode functions is K = 4.

K Value
The Decision Accuracy of Each IMF Component εK,i  From Figure 4, it can be seen that the simulated signal contains three main center frequencies: 170 Hz, 540 Hz, and 910 Hz. When the rotating machinery signal is weak or covered by strong background noise, the signal characteristics (impact property and center frequency) may not appear apparently. Comparing Figure 4 with Figure 5, it can be seen that the periodic impact components in the simulated signal are completely submerged by noise, the center frequency 910 Hz is difficult to distinguish, and the highfrequency part is severely affected by noise, with many interference spectral lines appear in this part.
In order to denoise the simulated signal, the VMD-MODWPT method is used. Firstly, the simulated signal z(t) is decomposed through VMD. The center frequency of each IMF component of z(t) is obtained under different K values. In order to obtain the optimal value of K, the judgment accuracy of each IMF component corresponding to different K values is calculated according to the method proposed in this paper. The calculation results are shown in Table 2. Table 2 shows that N4 = 4 > N5 = 1, so K = 4 can be tentatively set as the optimal value and N6 = N5 = 1. Therefore, the optimal number of mode functions is K = 4.

K Value
The Decision Accuracy of Each IMF Component εK,i From Figure 4, it can be seen that the simulated signal contains three main center frequencies: 170 Hz, 540 Hz, and 910 Hz. When the rotating machinery signal is weak or covered by strong background noise, the signal characteristics (impact property and center frequency) may not appear apparently. Comparing Figure 4 with Figure 5, it can be seen that the periodic impact components in the simulated signal are completely submerged by noise, the center frequency 910 Hz is difficult to distinguish, and the high-frequency part is severely affected by noise, with many interference spectral lines appear in this part.
In order to denoise the simulated signal, the VMD-MODWPT method is used. Firstly, the simulated signal z(t) is decomposed through VMD. The center frequency of each IMF component of z(t) is obtained under different K values. In order to obtain the optimal value of K, the judgment accuracy of each IMF component corresponding to different K values is calculated according to the method proposed in this paper. The calculation results are shown in Table 2. Table 2 shows that N 4 = 4 > N 5 = 1, so K = 4 can be tentatively set as the optimal value and N 6 = N 5 = 1. Therefore, the optimal number of mode functions is K = 4.

K Value
The Decision Accuracy of Each IMF Component ε K,i This method is applied to the rolling bearing vibration data of Case Western University in [40,41]. The analysis proves the effectiveness of the algorithm for selecting the value of K based on the number N of the effective center frequencies. Due to space limitations, only Table 1 is listed in [40]. The decision accuracies of each IMF component of the bearing inner ring fault signal corresponding to different K values are provided in Table 3. Table 3. The decision accuracy of each IMF component of the bearing inner ring fault signal.

K Value
The Decision Accuracy of Each IMF Component ε K,i From Table 3, it can be seen that N 2 = 1 > N 3 = 0, K = 2 can be tentatively set to the optimal value. At the same time, N 4 = 4 > N 5 = 3 and N 4 > N 2 . Therefore, the optimal value of the preset scale is K = 4. The above analysis verifies the effectiveness of the proposed method; it also proves the importance of step 8 in the algorithm for the signals composed of multiple frequency components. And if θ > 0.15 (for example θ = 0.16), this may cause N 5 = 6 > N 4 = 4, resulting in excessive selection of K and causing over-decomposition and modal mixing problems. Further, the analysis of the measured rolling bearing and the gearbox fault signals from [42,43] proves the effectiveness and practicality of this method.
When K = 4, the penalty factor α is calculated in different ranges, and the mean MDE value of the reconstruction signal is shown in Figure 6. During MDE calculation, the dimension of the embedded vector m = 3, the number of categories c = 6, the time delay d = 1, and the scale factor τ max = 10.

Sensors 2023, 23, x FOR PEER REVIEW
This method is applied to the rolling bearing vibration data of Case Western sity in [40,41]. The analysis proves the effectiveness of the algorithm for selecting th of K based on the number N of the effective center frequencies. Due to space limi only Table 1 is listed in [40]. The decision accuracies of each IMF component of the inner ring fault signal corresponding to different K values are provided in Table 3 Table 3, it can be seen that N2 = 1 > N3 = 0, K = 2 can be tentatively se optimal value. At the same time, N4 = 4 > N5 = 3 and N4 > N2. Therefore, the optim of the preset scale is K = 4. The above analysis verifies the effectiveness of the pr method; it also proves the importance of step 8 in the algorithm for the signals com of multiple frequency components. And if θ > 0.15 (for example θ = 0.16), this ma N5 = 6 > N4 = 4, resulting in excessive selection of K and causing over-decomposit modal mixing problems. Further, the analysis of the measured rolling bearing gearbox fault signals from [42,43] proves the effectiveness and practicality of this m When K = 4, the penalty factor α is calculated in different ranges, and the mea value of the reconstruction signal is shown in Figure 6. During MDE calculation mension of the embedded vector m = 3, the number of categories c = 6, the time d 1, and the scale factor τmax = 10. In Figure 6, when α = 3600, the mean MDE value of the reconstructed sign VMD decomposition is the smallest. This indicates that the impact component re the fault feature in the reconstructed signal contains the most, and also shows stro ularity and self-similarity, so α = 3600 is taken to decompose the simulation signa According to Table 1 and Figure 4, K = 4 and α = 3600 are used in VMD. Th decomposition results of the simulated signal z(t) and the spectrum of each IMF In Figure 6, when α = 3600, the mean MDE value of the reconstructed signal after VMD decomposition is the smallest. This indicates that the impact component related to the fault feature in the reconstructed signal contains the most, and also shows strong regularity and self-similarity, so α = 3600 is taken to decompose the simulation signal.
According to Table 1 and Figure 4, K = 4 and α = 3600 are used in VMD. The VMD decomposition results of the simulated signal z(t) and the spectrum of each IMF component are shown in Figure 7. To compare the decomposition effect, the EMD method [12] and the EEMD method are used to decompose the simulated signal z(t). The decomposition results and the spectrum of each IMF component are shown in Figures 8 and 9. An ensemble size of I = 100 and standard deviation ε0 = 0.2 are used in EEMD [13,44]. Further, for the sake of convenience or convenience of contrast, EMD and EEMD take the first to fourth IMF for analysis. It can be seen from Figure 7 that the VMD decomposition results are reasonable, in which the IMF1-IMF3 are amplitude modulation-frequency modulation signals, whereas the IMF4 is a high-frequency Gaussian white noise component. The frequency of the IMF1-IMF3 components is mainly concentrated near the center frequency. It is verified that the parameter selection method can effectively suppress the modal aliasing problem generated in the decomposition process. At the same time, it can reduce the information leakage between the modal components. However, it can also be seen that the frequency characteristics of the IMF3 component are hidden by noise, so its feature information must be effectively extracted and restored.
To compare the decomposition effect, the EMD method [12] and the EEMD method are used to decompose the simulated signal z(t). The decomposition results and the spectrum of each IMF component are shown in Figures 8 and 9. An ensemble size of I = 100 and standard deviation ε 0 = 0.2 are used in EEMD [13,44]. Further, for the sake of convenience or convenience of contrast, EMD and EEMD take the first to fourth IMF for analysis. To compare the decomposition effect, the EMD method [12] and the EEMD method are used to decompose the simulated signal z(t). The decomposition results and the spectrum of each IMF component are shown in Figures 8 and 9. An ensemble size of I = 100 and standard deviation ε0 = 0.2 are used in EEMD [13,44]. Further, for the sake of convenience or convenience of contrast, EMD and EEMD take the first to fourth IMF for analysis.    According to the information given in Figures 6 and 7, EEMD suppresses the modal mixing problem to a certain extent. However, it can also be seen that the decomposition result is not ideal. EMD and EEMD obtain many iteration error components, especially for the IMF1, in which almost all the frequency components exist in its entire frequency band. This inevitably affects the accuracy of the subsequent signal feature extractions.
The CWE value of each IMF component decomposed through VMD is calculated and the principal component factor λ i is obtained based on CWE. The λ i of each IMF component is show in Figure 10. During CWE calculation, δ = 0.2 and β = γ = 0.4. s 2023, 23, x FOR PEER REVIEW According to the information given in Figures 6 and 7, EEMD mixing problem to a certain extent. However, it can also be seen result is not ideal. EMD and EEMD obtain many iteration error for the IMF1, in which almost all the frequency components exis band. This inevitably affects the accuracy of the subsequent signa The CWE value of each IMF component decomposed through the principal component factor λi is obtained based on CWE. The nent is show in Figure 10. During CWE calculation, δ = 0.2 and β =  Figure 10, it can be seen that the maximum difference nent factors λi is between IMF1 and IMF3, so the first two compo after reordering are the components, which containing the main c nal. Therefore, the illusive IMF component selection algorith properly remove the interference components. According to the p is a high-frequency amplitude modulation-frequency modulation ulation signal, whose characteristics are hidden by noise. In order effect, it is necessary to extract its feature. For the IMF4 compone nent factor λ4 = 0, so this component is recognized as the noise com eliminated.
In the decomposition process, the IMF3 component is MOWDPT. Based on [28], during the decomposition process, the From Figure 10, it can be seen that the maximum difference of the principal component factors λ i is between IMF1 and IMF3, so the first two components (IMF2 and IMF1) after reordering are the components, which containing the main characteristics of the signal. Therefore, the illusive IMF component selection algorithm based on CWE can properly remove the interference components. According to the previous analysis, IMF3 is a highfrequency amplitude modulation-frequency modulation component of the simulation signal, whose characteristics are hidden by noise. In order to ensure the denoising effect, it is necessary to extract its feature. For the IMF4 component, its principal component factor λ 4 = 0, so this component is recognized as the noise component, which shall be eliminated.
In the decomposition process, the IMF3 component is decomposed through MOWDPT. Based on [28], during the decomposition process, the Fejer-Korovkin wavelet filter with length L = 22 and decomposition layer number J = 2 is selected. The decomposition results are shown in Figure 11.

eliminated.
In the decomposition process, the IMF3 component is decomposed thro MOWDPT. Based on [28], during the decomposition process, the Fejer-Korovkin wav filter with length L = 22 and decomposition layer number J = 2 is selected. The decom sition results are shown in Figure 11. In Figure 11, the impact characteristics of the C2 and C3 components are obvious the peak value of the spectral line at the central frequency of the C2 and C3 compon is large. In order to determine the component containing the characteristic informatio In Figure 11, the impact characteristics of the C2 and C3 components are obvious and the peak value of the spectral line at the central frequency of the C2 and C3 components is large. In order to determine the component containing the characteristic information of the signal, the value of the DID between each component and the normal signal is calculated. The results are shown in Table 4. In DID calculation, α 1 = 100, α 2 = 0.1. It can be seen from Table 3 that the value of the DID between the C2 component and the normal signal is the smallest. This indicates that the above components contain lots of information on the normal characteristics. Therefore, this component is judged as the high-frequency denoising component, and the denoising simulation signal is reconstructed using the IMF1 and IMF2 obtained through VMD and the C2 component of IMF3 obtained through MODWPT. The time domain waveform and the spectrum of the denoising signal are shown in Figure 12.  Table 4. In DID calculation, α1 = 100, α2 = 0.1. It can be seen from Table 3 that the value of the DID between the C2 component an the normal signal is the smallest. This indicates that the above components contain lots information on the normal characteristics. Therefore, this component is judged as t high-frequency denoising component, and the denoising simulation signal is reco structed using the IMF1 and IMF2 obtained through VMD and the C2 component of IM obtained through MODWPT. The time domain waveform and the spectrum of the d noising signal are shown in Figure 12. Comparing Figure 4 with Figure 12, it can be seen that after processing with the pr posed VMD-MODWPT method, most of the useless interference and noise componen in the simulated signal z(t) are filtered out effectively, and the signal time domain wavefor and the spectrum effectively highlight the characteristic information of the signal x(t).
To highlight the effectiveness and superiority of the proposed method, WTD, EM Comparing Figure 4 with Figure 12, it can be seen that after processing with the proposed VMD-MODWPT method, most of the useless interference and noise components in the simulated signal z(t) are filtered out effectively, and the signal time domain waveform and the spectrum effectively highlight the characteristic information of the signal x(t).
To highlight the effectiveness and superiority of the proposed method, WTD, EMD-RD [45], and EEMD-WTD [37] were used to deal with the simulation signal z(t). The denoising results of different methods are shown in Figures 13-15. Based on [13,44,46], WTD adopts the 'db3' wavelet basis function for three-layer decomposition and the 'Rigrsure' rule for soft threshold denoising; for EEMD-WTD, an ensemble size of I = 100 and standard deviation ε 0 = 0.2 are used in EEMD. Comparing Figure 4 with Figure 12, it can be seen that after processing with the pr posed VMD-MODWPT method, most of the useless interference and noise componen in the simulated signal z(t) are filtered out effectively, and the signal time domain wavefor and the spectrum effectively highlight the characteristic information of the signal x(t).
To highlight the effectiveness and superiority of the proposed method, WTD, EM RD [45], and EEMD-WTD [37] were used to deal with the simulation signal z(t). The d noising results of different methods are shown in Figures 13-15. Based on [13,44,46], WT adopts the 'db3′ wavelet basis function for three-layer decomposition and the 'Rigrsu rule for soft threshold denoising; for EEMD-WTD, an ensemble size of I = 100 and stan ard deviation ε0 = 0.2 are used in EEMD. It can be seen from Figures 13-15 that the denoising effect of the above method is n ideal, because a large amount of noise remains in the signal even after denoising, whi inevitably affects the subsequent signal analysis effect. For the WTD method, the ma characteristic information of the simulated signal is also filtered. The center frequenc 540 Hz and 910 Hz do not show spectral lines in the spectrum. For the EMD-RD metho there is no obvious spectral line of the center frequency 170 Hz in the spectrum, and ma interference noises appear in the high-frequency part. Although the center frequencies the simulated signal are obtained through the EEMD-WTD method, lots of noise still a peared in the spectrum and time domain waveform of the denoising signal.
To quantitatively evaluate the noise reduction performance of the denoising me ods, the root mean squared error (RMSE) and the peak signal-to-noise ratio (PSNR) a It can be seen from Figures 13-15 that the denoising effect of the above method is n ideal, because a large amount of noise remains in the signal even after denoising, wh inevitably affects the subsequent signal analysis effect. For the WTD method, the m characteristic information of the simulated signal is also filtered. The center frequenc 540 Hz and 910 Hz do not show spectral lines in the spectrum. For the EMD-RD meth there is no obvious spectral line of the center frequency 170 Hz in the spectrum, and ma interference noises appear in the high-frequency part. Although the center frequencies the simulated signal are obtained through the EEMD-WTD method, lots of noise still a peared in the spectrum and time domain waveform of the denoising signal.
To quantitatively evaluate the noise reduction performance of the denoising me ods, the root mean squared error (RMSE) and the peak signal-to-noise ratio (PSNR) It can be seen from Figures 13-15 that the denoising effect of the above method is not ideal, because a large amount of noise remains in the signal even after denoising, which inevitably affects the subsequent signal analysis effect. For the WTD method, the main characteristic information of the simulated signal is also filtered. The center frequencies 540 Hz and 910 Hz do not show spectral lines in the spectrum. For the EMD-RD method, there is no obvious spectral line of the center frequency 170 Hz in the spectrum, and many interference noises appear in the high-frequency part. Although the center frequencies of the simulated signal are obtained through the EEMD-WTD method, lots of noise still appeared in the spectrum and time domain waveform of the denoising signal.
To quantitatively evaluate the noise reduction performance of the denoising methods, the root mean squared error (RMSE) and the peak signal-to-noise ratio (PSNR) are used as the evaluation factors. The corresponding calculation formulae of PSNR and PMSE are described as follows: where x(t) is the original signal, x (t) is the denoised signal, and N is the length of the signal. The calculation results of the evaluation factors of the denoising effect of different methods are shown in Table 5. Compared with WTD, EMD-RD, and EEMD-WTD, the evaluation factors of the denoising signal by using VMD-MODWPT are better than other three methods. Combined with the analysis in Figures 12-15, the effectiveness of the proposed method is proven.
The denoising evaluation factors of RMSE and PSNR from the above four methods for the simulated signal z(t) with different SNRs from −20 dB to 10 dB obtained by using the 'awgn' function in MATLAB software are shown in Table 6. As shown in Table 6, the proposed VMD-MODWPT method can denoise the signals effectively. Compared with WTD, EMD-RD, and EEMD-WTD, VMD-MODWPT shows the best performance of RMSE and PSNR at different SNRs.

Experimental Case Study
Gears are installed in different kinds of machinery. Several problems of these machines may be caused due to defects within them. Generally, in the fault diagnosis of rotating machines, the signals to be analyzed often come with a lot of noise. To show the efficiency of the proposed method, in this section, the gear fault signal, which is mainly used to simulate the related faults of common mechanical equipment, is selected for analysis.
The test apparatus used in this study is shown in Figure 16. The type of gear tester used is a QPZZ-II and the experiment set-up consists of a single-stage gearbox driven by a 0.75 kW AC governor motor. The driving gear consists of 55 teeth and the driven gear of 75 teeth, respectively. The signal monitoring and collection system consists of a KD1001L accelerometer, an amplifier, a data collector, and a computer. The sampling frequency is 5120 Hz. After sampling, the collected vibration signals are loaded into MATLAB from the data files. During the test, the driving motor's rotating frequency is 50 Hz. The acce sensor is installed on the bearing seat of the gearbox. The tape of the data acquisit is ADA16-8/2 (LPCI). Each vibration signal consists of 5120 data points and the s precision is 16 bt. The parameters of the gearbox are shown in Table 7. The pitting fault is the cause of the fault of the driving gear. Figure 17 sh driving gear; the main cause of this fault is the spot welding on the gear tooth su makes the welding slag peel off, then a pit is formed on the tooth surface. During the test, the driving motor's rotating frequency is 50 Hz. The acceleration sensor is installed on the bearing seat of the gearbox. The tape of the data acquisition card is ADA16-8/2 (LPCI). Each vibration signal consists of 5120 data points and the sampling precision is 16 bt. The parameters of the gearbox are shown in Table 7. The pitting fault is the cause of the fault of the driving gear. Figure 17 shows the driving gear; the main cause of this fault is the spot welding on the gear tooth surface. It makes the welding slag peel off, then a pit is formed on the tooth surface. The pitting fault is the cause of the fault of the driving gear. Figure 17 sh driving gear; the main cause of this fault is the spot welding on the gear tooth su makes the welding slag peel off, then a pit is formed on the tooth surface. For gear fault diagnosis, envelope spectrum analysis is the most common method. Therefore, it is necessary to resolve the envelope spectrum of the gear fau to diagnose the type of gear fault. The time domain waveform and envelope spec the gear fault signal are shown in Figure 18. For gear fault diagnosis, envelope spectrum analysis is the most commonly used method. Therefore, it is necessary to resolve the envelope spectrum of the gear fault signal to diagnose the type of gear fault. The time domain waveform and envelope spectrum of the gear fault signal are shown in Figure 18. It can be seen in Figure 18a,b that because no corresponding noise cancellation dev is used in the signal acquisition system, the collected time domain signal contains mo noise components. The maximum envelope spectrum is focused in 0-300 Hz of the f quency band. For ease of analysis, the envelope spectrum limited to 250 Hz is shown Figure 16c. However, the fault characteristics and fault types of gears can be read a distinguished with difficulty; the main reason for the above problems is the influence the noise interference components.
The time domain waveform and envelope spectrum of the denoised signal using t proposed VMD-MODWPT method are shown in Figure 19. It can be seen in Figure 18a,b that because no corresponding noise cancellation device is used in the signal acquisition system, the collected time domain signal contains more noise components. The maximum envelope spectrum is focused in 0-300 Hz of the frequency band. For ease of analysis, the envelope spectrum limited to 250 Hz is shown in Figure 16c. However, the fault characteristics and fault types of gears can be read and distinguished with difficulty; the main reason for the above problems is the influence of the noise interference components.
The time domain waveform and envelope spectrum of the denoised signal using the proposed VMD-MODWPT method are shown in Figure 19. quency band. For ease of analysis, the envelope spectrum limited to 250 Hz is shown in Figure 16c. However, the fault characteristics and fault types of gears can be read and distinguished with difficulty; the main reason for the above problems is the influence of the noise interference components.
The time domain waveform and envelope spectrum of the denoised signal using the proposed VMD-MODWPT method are shown in Figure 19. It can be seen from Figure 19 that most of the useless interference and noise components are filtered out effectively. In the envelope spectrum of the gear pitting fault denoising signal, 14.67 Hz and 29.38 Hz have obvious spectral lines, and 14.67 Hz and 29.38 Hz are about one and two times the rotating frequency of the driving gear. According to the principle of gear vibration [47], it is known that when a gear has a pitting fault, the rotating frequency and frequency doubling components of the shaft where the faulty gear is located can be obtained through envelope spectrum analysis, so it can be determined that the driving gear has a pitting fault.
To highlight the effectiveness and superiority of the proposed method, the WTD method, EMD-RD method, and EEMD-WTD method are used to deal with the gear pitting fault signal. The denoising results of the different methods are shown in Figures 20-22. It can be seen from Figure 19 that most of the useless interference and noise components are filtered out effectively. In the envelope spectrum of the gear pitting fault denoising signal, 14.67 Hz and 29.38 Hz have obvious spectral lines, and 14.67 Hz and 29.38 Hz are about one and two times the rotating frequency of the driving gear. According to the principle of gear vibration [47], it is known that when a gear has a pitting fault, the rotating frequency and frequency doubling components of the shaft where the faulty gear is located can be obtained through envelope spectrum analysis, so it can be determined that the driving gear has a pitting fault.
To highlight the effectiveness and superiority of the proposed method, the WTD method, EMD-RD method, and EEMD-WTD method are used to deal with the gear pitting fault signal. The denoising results of the different methods are shown in . Based on [13,43,44,46], WTD adopts the 'db3' wavelet basis function for three-layer decomposition and the 'Rigrsure' rule for soft threshold denoising; for EEMD-WTD, an ensemble size of I = 100 and standard deviation ε 0 = 0.2 are used in EEMD. Based on [13,43,44,46], WTD adopts the 'db3′ wavelet basis function for three-layer decomposition and the 'Rigrsure' rule for soft threshold denoising; for EEMD-WTD, an ensemble size of I = 100 and standard deviation ε0 = 0.2 are used in EEMD.  Based on [13,43,44,46], WTD adopts the 'db3′ wavelet basis function for three-layer decomposition and the 'Rigrsure' rule for soft threshold denoising; for EEMD-WTD, an ensemble size of I = 100 and standard deviation ε0 = 0.2 are used in EEMD.   In Figures 20-22, due to the influence of the noise and mode mixing problem, there are many interfering frequency components in the envelope spectra of the denoising signal obtained using WTD, EMD-RD, and EEMD-WTD. The spectral line features of 14.67 Hz and 29.38 Hz are not obvious, but 10.63 Hz and 21.25 Hz are clearly observed, which would reduce the accuracy of the fault feature extraction and the reliability of fault diagnosis, with 10.63 Hz being close to the rotating frequency of the driven gear (10.76 Hz). It also means that the completeness and denoising effect of the proposed method are much better than the other three decomposition methods.
In order to verify the capability of the VMD-MODWPT method for dealing with the rotating machinery signals with strong noise, a Gaussian noise signal with −15 dB of SNR is added to the gear pitting fault signal, which is obtained by using the 'awgn' function in MATLAB software. Its time domain waveform and envelope spectrum are shown in  . It also means that the completeness and denoising effect of the proposed method are much better than the other three decomposition methods.
In order to verify the capability of the VMD-MODWPT method for dealing with the rotating machinery signals with strong noise, a Gaussian noise signal with −15 dB of SNR is added to the gear pitting fault signal, which is obtained by using the 'awgn' function in MATLAB software. Its time domain waveform and envelope spectrum are shown in Figure 23. It can be seen from Figure 23 that the impulsive components in the time doma waveform are submerged by the noise. The envelope spectrum in the whole frequen band does not have obvious characteristics. The maximum envelope spectrum is focus in the 0-200 Hz frequency band. However, in this frequency band, it is also difficult identify the fault characteristics.
The VMD-MODWPT method is used to denoise the gear pitting fault signal w strong noise. The time domain waveform and envelope spectrum of the denoised sign are shown in Figure 22. It can be seen from Figure 23 that the impulsive components in the time domain waveform are submerged by the noise. The envelope spectrum in the whole frequency band does not have obvious characteristics. The maximum envelope spectrum is focused in the 0-200 Hz frequency band. However, in this frequency band, it is also difficult to identify the fault characteristics.
The VMD-MODWPT method is used to denoise the gear pitting fault signal with strong noise. The time domain waveform and envelope spectrum of the denoised signal are shown in Figure 22.
Comparing the time domain waveform and envelope spectrum of the denoised signal in Figure 24 with those of Figure 23, it can be seen that the impulsive components are obvious in the denoised signal, and 14.38 Hz, 29.38 Hz, an d76.25 Hz have obvious spectral lines in the envelope spectrum, which are the rotating frequency of the driving gear and multiples of its frequency. Therefore, under the condition of strong noise, it is able to effectively diagnose the pitting fault of gears.
waveform are submerged by the noise. The envelope spectrum in the whole frequen band does not have obvious characteristics. The maximum envelope spectrum is focus in the 0-200 Hz frequency band. However, in this frequency band, it is also difficult identify the fault characteristics.
The VMD-MODWPT method is used to denoise the gear pitting fault signal w strong noise. The time domain waveform and envelope spectrum of the denoised sig are shown in Figure 22.
Comparing the time domain waveform and envelope spectrum of the denoised s nal in Figure 24 with those of Figure 23, it can be seen that the impulsive components obvious in the denoised signal, and 14.38 Hz, 29.38 Hz, an d76.25 Hz have obvious spect lines in the envelope spectrum, which are the rotating frequency of the driving gear a multiples of its frequency. Therefore, under the condition of strong noise, it is able to fectively diagnose the pitting fault of gears. The gear pitting fault signal with strong noise is also decomposed with the WT method, the EMD-RD method, and the EEMD-WTD method. The denoised results shown in Figures 25-27. Based on [13,[44][45][46], WTD adopts the 'db3′ wavelet basis functi for three-layer decomposition and the 'Rigrsure' rule for soft threshold denoising; EEMD-WTD, an ensemble size of I = 100 and standard deviation ε0 = 0.2 are used in EEM The gear pitting fault signal with strong noise is also decomposed with the WTD method, the EMD-RD method, and the EEMD-WTD method. The denoised results are shown in Figures 25-27. Based on [13,[44][45][46], WTD adopts the 'db3' wavelet basis function for three-layer decomposition and the 'Rigrsure' rule for soft threshold denoising; for EEMD-WTD, an ensemble size of I = 100 and standard deviation ε 0 = 0.2 are used in EEMD.    As can be seen from Figures 25-27, the noise in the signal has not been effectiv filtered, which means some impulsive components are hidden by the noise. Because lots of interfering frequency components in the envelope spectrum of the denoised sign the spectral line features of 14.67 Hz and 29.38 Hz are not obvious. The results proved th the completeness and denoising effect of the proposed method are much better than t other three decomposition methods, even for the signal with strong noise.
In order to quantify the evaluation of the denoising effect, the RMSE and PSNR calculated for the four methods. The results are shown in Table 8 (for convenience of c culation, the collected gear pitting fault signal is treated as the clear signal). The resu show that even for the signal with strong noise, the proposed VMD-MODWPT meth can denoise with better performance.  As can be seen from Figures 25-27 the noise in the signal has not been effectively filtered, which means some impulsive components are hidden by the noise. Because of lots of interfering frequency components in the envelope spectrum of the denoised signal, the spectral line features of 14.67 Hz and 29.38 Hz are not obvious. The results proved that the completeness and denoising effect of the proposed method are much better than the other three decomposition methods, even for the signal with strong noise.
In order to quantify the evaluation of the denoising effect, the RMSE and PSNR are calculated for the four methods. The results are shown in Table 8 (for convenience of calculation, the collected gear pitting fault signal is treated as the clear signal). The results show that even for the signal with strong noise, the proposed VMD-MODWPT method can denoise with better performance.

Conclusions
This paper proposes a signal denoising method based on a fault feature of rotating machinery. The principle of minimum MDE and the number of effective centers N are proposed as parameters to be used in the selection of parameters for the signal. The CWE value is used to calculate the components that represent decomposition errors or background signals. The components with high frequencies are denoised and reconstructed. The denoised signal is reconstructed with sensitive and high-frequency components. The proposed method has advanced denoising performance when analyzing both simulation and real gear pitting fault signals, regardless of whether the signal has a high or low signal-to-noise ratio.
In this paper, the rotating machinery signal is investigated. It is thought that the method described in this paper performs better for other signals. In the future, we might examine the power system signal using this approach and evaluate it against other signal denoising methods. In addition, we will also try to study whether the proposed method is suitable only for steady states or whether it can also handle transitions or dynamic behaviour of rotary machinery.